Usenet 1993 July/InfoMagic USENET CD-ROM July 1993.ISO/sources/misc/volume24/gnucalc/part10

Text File | 1991-10-29 | 55.3 KB | 1,789 lines

Newsgroups: comp.sources.misc
From: [e-mail] (David Gillespie)
Subject:  v24i058:  gnucalc - GNU Emacs Calculator, v2.00, Part10/56
Message-ID: <[e-mail]>
X-Md4-Signature: 017f2514de54b02b0b894f816e81dee4
Date: Tue, 29 Oct 1991 23:00:26 GMT
Approved: [e-mail]

Submitted-by: [e-mail] (David Gillespie)
Posting-number: Volume 24, Issue 58
Archive-name: gnucalc/part10
Environment: Emacs
Supersedes: gmcalc: Volume 13, Issue 27-45

---- Cut Here and unpack ----
#!/bin/sh
# this is Part.10 (part 10 of a multipart archive)
# do not concatenate these parts, unpack them in order with /bin/sh
# file calc-bin.el continued
#
if test ! -r _shar_seq_.tmp; then
    echo 'Please unpack part 1 first!'
    exit 1
fi
(read Scheck
 if test "$Scheck" != 10; then
    echo Please unpack part "$Scheck" next!
    exit 1
 else
    exit 0
 fi
) < _shar_seq_.tmp || exit 1
if test ! -f _shar_wnt_.tmp; then
    echo 'x - still skipping calc-bin.el'
else
echo 'x - continuing file calc-bin.el'
sed 's/^X//' << 'SHAR_EOF' >> 'calc-bin.el' &&
X  (interactive "NDisplay radix (2-36): ")
X  (calc-wrapper
X   (if (and (>= n 2) (<= n 36))
X       (progn
X     (calc-change-mode 'calc-number-radix n t)
X     ;; also change global value so minibuffer sees it
X     (setq-default calc-number-radix calc-number-radix))
X     (setq n calc-number-radix))
X   (message "Number radix is %d." n))
)
X
(defun calc-decimal-radix ()
X  (interactive)
X  (calc-radix 10)
)
X
(defun calc-binary-radix ()
X  (interactive)
X  (calc-radix 2)
)
X
(defun calc-octal-radix ()
X  (interactive)
X  (calc-radix 8)
)
X
(defun calc-hex-radix ()
X  (interactive)
X  (calc-radix 16)
)
X
(defun calc-leading-zeros (n)
X  (interactive "P")
X  (calc-wrapper
X   (if (calc-change-mode 'calc-leading-zeros n t t)
X       (message "Zero-padding integers to %d digits (assuming radix %d)."
X        (let* ((calc-internal-prec 6))
X          (math-compute-max-digits (math-abs calc-word-size)
X                       calc-number-radix))
X        calc-number-radix)
X     (message "Omitting leading zeros on integers.")))
)
X
X
(defvar math-power-of-2-cache (list 1 2 4 8 16 32 64 128 256 512 1024))
(defvar math-big-power-of-2-cache nil)
(defun math-power-of-2 (n)    ;  [I I] [Public]
X  (if (and (natnump n) (<= n 100))
X      (or (nth n math-power-of-2-cache)
X      (let* ((i (length math-power-of-2-cache))
X         (val (nth (1- i) math-power-of-2-cache)))
X        (while (<= i n)
X          (setq val (math-mul val 2)
X            math-power-of-2-cache (nconc math-power-of-2-cache
X                         (list val))
X            i (1+ i)))
X        val))
X    (let ((found (assq n math-big-power-of-2-cache)))
X      (if found
X      (cdr found)
X    (let ((po2 (math-ipow 2 n)))
X      (setq math-big-power-of-2-cache
X        (cons (cons n po2) math-big-power-of-2-cache))
X      po2))))
)
X
(defun math-integer-log2 (n)    ; [I I] [Public]
X  (let ((i 0)
X    (p math-power-of-2-cache)
X    val)
X    (while (and p (Math-natnum-lessp (setq val (car p)) n))
X      (setq p (cdr p)
X        i (1+ i)))
X    (if p
X    (and (equal val n)
X         i)
X      (while (Math-natnum-lessp
X          (prog1
X          (setq val (math-mul val 2))
X        (setq math-power-of-2-cache (nconc math-power-of-2-cache
X                           (list val))))
X          n)
X    (setq i (1+ i)))
X      (and (equal val n)
X       i)))
)
X
X
X
X
;;; Bitwise operations.
X
(defun calcFunc-and (a b &optional w)   ; [I I I] [Public]
X  (cond ((Math-messy-integerp w)
X     (calcFunc-and a b (math-trunc w)))
X    ((and w (not (integerp w)))
X     (math-reject-arg w 'fixnump))
X    ((and (integerp a) (integerp b))
X     (math-clip (logand a b) w))
X    ((or (eq (car-safe a) 'mod) (eq (car-safe b) 'mod))
X     (math-binary-modulo-args 'calcFunc-and a b w))
X    ((not (Math-num-integerp a))
X     (math-reject-arg a 'integerp))
X    ((not (Math-num-integerp b))
X     (math-reject-arg b 'integerp))
X    (t (math-clip (cons 'bigpos
X                (math-and-bignum (math-binary-arg a w)
X                         (math-binary-arg b w)))
X              w)))
)
X
(defun math-binary-arg (a w)
X  (if (not (Math-integerp a))
X      (setq a (math-trunc a)))
X  (if (Math-integer-negp a)
X      (math-not-bignum (cdr (math-bignum-test (math-sub -1 a)))
X               (math-abs (if w (math-trunc w) calc-word-size)))
X    (cdr (Math-bignum-test a)))
)
X
(defun math-binary-modulo-args (f a b w)
X  (let (mod)
X    (if (eq (car-safe a) 'mod)
X    (progn
X      (setq mod (nth 2 a)
X        a (nth 1 a))
X      (if (eq (car-safe b) 'mod)
X          (if (equal mod (nth 2 b))
X          (setq b (nth 1 b))
X        (math-reject-arg b "*Inconsistent modulos"))))
X      (setq mod (nth 2 b)
X        b (nth 1 b)))
X    (if (Math-messy-integerp mod)
X    (setq mod (math-trunc mod))
X      (or (Math-integerp mod)
X      (math-reject-arg mod 'integerp)))
X    (let ((bits (math-integer-log2 mod)))
X      (if bits
X      (if w
X          (if (/= w bits)
X          (calc-record-why
X           "*Warning: Modulo inconsistent with word size"))
X        (setq w bits))
X    (calc-record-why "*Warning: Modulo is not a power of 2"))
X      (math-make-mod (if b
X             (funcall f a b w)
X               (funcall f a w))
X             mod)))
)
X
(defun math-and-bignum (a b)   ; [l l l]
X  (and a b
X       (let ((qa (math-div-bignum-digit a 512))
X         (qb (math-div-bignum-digit b 512)))
X     (math-mul-bignum-digit (math-and-bignum (math-norm-bignum (car qa))
X                          (math-norm-bignum (car qb)))
X                 512
X                 (logand (cdr qa) (cdr qb)))))
)
X
(defun calcFunc-or (a b &optional w)   ; [I I I] [Public]
X  (cond ((Math-messy-integerp w)
X     (calcFunc-or a b (math-trunc w)))
X    ((and w (not (integerp w)))
X     (math-reject-arg w 'fixnump))
X    ((and (integerp a) (integerp b))
X     (math-clip (logior a b) w))
X    ((or (eq (car-safe a) 'mod) (eq (car-safe b) 'mod))
X     (math-binary-modulo-args 'calcFunc-or a b w))
X    ((not (Math-num-integerp a))
X     (math-reject-arg a 'integerp))
X    ((not (Math-num-integerp b))
X     (math-reject-arg b 'integerp))
X    (t (math-clip (cons 'bigpos
X                (math-or-bignum (math-binary-arg a w)
X                        (math-binary-arg b w)))
X              w)))
)
X
(defun math-or-bignum (a b)   ; [l l l]
X  (and (or a b)
X       (let ((qa (math-div-bignum-digit a 512))
X         (qb (math-div-bignum-digit b 512)))
X     (math-mul-bignum-digit (math-or-bignum (math-norm-bignum (car qa))
X                         (math-norm-bignum (car qb)))
X                 512
X                 (logior (cdr qa) (cdr qb)))))
)
X
(defun calcFunc-xor (a b &optional w)   ; [I I I] [Public]
X  (cond ((Math-messy-integerp w)
X     (calcFunc-xor a b (math-trunc w)))
X    ((and w (not (integerp w)))
X     (math-reject-arg w 'fixnump))
X    ((and (integerp a) (integerp b))
X     (math-clip (logxor a b) w))
X    ((or (eq (car-safe a) 'mod) (eq (car-safe b) 'mod))
X     (math-binary-modulo-args 'calcFunc-xor a b w))
X    ((not (Math-num-integerp a))
X     (math-reject-arg a 'integerp))
X    ((not (Math-num-integerp b))
X     (math-reject-arg b 'integerp))
X    (t (math-clip (cons 'bigpos
X                (math-xor-bignum (math-binary-arg a w)
X                         (math-binary-arg b w)))
X              w)))
)
X
(defun math-xor-bignum (a b)   ; [l l l]
X  (and (or a b)
X       (let ((qa (math-div-bignum-digit a 512))
X         (qb (math-div-bignum-digit b 512)))
X     (math-mul-bignum-digit (math-xor-bignum (math-norm-bignum (car qa))
X                          (math-norm-bignum (car qb)))
X                 512
X                 (logxor (cdr qa) (cdr qb)))))
)
X
(defun calcFunc-diff (a b &optional w)   ; [I I I] [Public]
X  (cond ((Math-messy-integerp w)
X     (calcFunc-diff a b (math-trunc w)))
X    ((and w (not (integerp w)))
X     (math-reject-arg w 'fixnump))
X    ((and (integerp a) (integerp b))
X     (math-clip (logand a (lognot b)) w))
X    ((or (eq (car-safe a) 'mod) (eq (car-safe b) 'mod))
X     (math-binary-modulo-args 'calcFunc-diff a b w))
X    ((not (Math-num-integerp a))
X     (math-reject-arg a 'integerp))
X    ((not (Math-num-integerp b))
X     (math-reject-arg b 'integerp))
X    (t (math-clip (cons 'bigpos
X                (math-diff-bignum (math-binary-arg a w)
X                          (math-binary-arg b w)))
X              w)))
)
X
(defun math-diff-bignum (a b)   ; [l l l]
X  (and a
X       (let ((qa (math-div-bignum-digit a 512))
X         (qb (math-div-bignum-digit b 512)))
X     (math-mul-bignum-digit (math-diff-bignum (math-norm-bignum (car qa))
X                           (math-norm-bignum (car qb)))
X                 512
X                 (logand (cdr qa) (lognot (cdr qb))))))
)
X
(defun calcFunc-not (a &optional w)   ; [I I] [Public]
X  (cond ((Math-messy-integerp w)
X     (calcFunc-not a (math-trunc w)))
X    ((eq (car-safe a) 'mod)
X     (math-binary-modulo-args 'calcFunc-not a nil w))
X    ((and w (not (integerp w)))
X     (math-reject-arg w 'fixnump))
X    ((not (Math-num-integerp a))
X     (math-reject-arg a 'integerp))
X    ((< (or w (setq w calc-word-size)) 0)
X     (math-clip (calcFunc-not a (- w)) w))
X    (t (math-normalize
X        (cons 'bigpos
X          (math-not-bignum (math-binary-arg a w)
X                   w)))))
)
X
(defun math-not-bignum (a w)   ; [l l]
X  (let ((q (math-div-bignum-digit a 512)))
X    (if (<= w 9)
X    (list (logand (lognot (cdr q))
X              (1- (lsh 1 w))))
X      (math-mul-bignum-digit (math-not-bignum (math-norm-bignum (car q))
X                           (- w 9))
X                  512
X                  (logxor (cdr q) 511))))
)
X
(defun calcFunc-lsh (a &optional n w)   ; [I I] [Public]
X  (setq a (math-trunc a)
X    n (if n (math-trunc n) 1))
X  (if (eq (car-safe a) 'mod)
X      (math-binary-modulo-args 'calcFunc-lsh a n w)
X    (setq w (if w (math-trunc w) calc-word-size))
X    (or (integerp w)
X    (math-reject-arg w 'fixnump))
X    (or (Math-integerp a)
X    (math-reject-arg a 'integerp))
X    (or (Math-integerp n)
X    (math-reject-arg n 'integerp))
X    (if (< w 0)
X    (math-clip (calcFunc-lsh a n (- w)) w)
X      (if (Math-integer-negp a)
X      (setq a (math-clip a w)))
X      (cond ((or (Math-lessp n (- w))
X         (Math-lessp w n))
X         0)
X        ((< n 0)
X         (math-quotient (math-clip a w) (math-power-of-2 (- n))))
X        (t
X         (math-clip (math-mul a (math-power-of-2 n)) w)))))
)
X
(defun calcFunc-rsh (a &optional n w)   ; [I I] [Public]
X  (calcFunc-lsh a (math-neg (or n 1)) w)
)
X
(defun calcFunc-ash (a &optional n w)   ; [I I] [Public]
X  (if (or (null n)
X      (not (Math-negp n)))
X      (calcFunc-lsh a n w)
X    (setq a (math-trunc a)
X      n (if n (math-trunc n) 1))
X    (if (eq (car-safe a) 'mod)
X    (math-binary-modulo-args 'calcFunc-ash a n w)
X      (setq w (if w (math-trunc w) calc-word-size))
X      (or (integerp w)
X      (math-reject-arg w 'fixnump))
X      (or (Math-integerp a)
X      (math-reject-arg a 'integerp))
X      (or (Math-integerp n)
X      (math-reject-arg n 'integerp))
X      (if (< w 0)
X      (math-clip (calcFunc-ash a n (- w)) w)
X    (if (Math-integer-negp a)
X        (setq a (math-clip a w)))
X    (let ((two-to-sizem1 (math-power-of-2 (1- w)))
X          (sh (calcFunc-lsh a n w)))
X      (cond ((Math-natnum-lessp a two-to-sizem1)
X         sh)
X        ((Math-lessp n (- 1 w))
X         (math-add (math-mul two-to-sizem1 2) -1))
X        (t (let ((two-to-n (math-power-of-2 (- n))))
X             (math-add (calcFunc-lsh (math-add two-to-n -1)
X                         (+ w n) w)
X                   sh))))))))
)
X
(defun calcFunc-rash (a &optional n w)   ; [I I] [Public]
X  (calcFunc-ash a (math-neg (or n 1)) w)
)
X
(defun calcFunc-rot (a &optional n w)   ; [I I] [Public]
X  (setq a (math-trunc a)
X    n (if n (math-trunc n) 1))
X  (if (eq (car-safe a) 'mod)
X      (math-binary-modulo-args 'calcFunc-rot a n w)
X    (setq w (if w (math-trunc w) calc-word-size))
X    (or (integerp w)
X    (math-reject-arg w 'fixnump))
X    (or (Math-integerp a)
X    (math-reject-arg a 'integerp))
X    (or (Math-integerp n)
X    (math-reject-arg n 'integerp))
X    (if (< w 0)
X    (math-clip (calcFunc-rot a n (- w)) w)
X      (if (Math-integer-negp a)
X      (setq a (math-clip a w)))
X      (cond ((or (Math-integer-negp n)
X         (not (Math-natnum-lessp n w)))
X         (calcFunc-rot a (math-mod n w) w))
X        (t
X         (math-add (calcFunc-lsh a (- n w) w)
X               (calcFunc-lsh a n w))))))
)
X
(defun math-clip (a &optional w)   ; [I I] [Public]
X  (cond ((Math-messy-integerp w)
X     (math-clip a (math-trunc w)))
X    ((eq (car-safe a) 'mod)
X     (math-binary-modulo-args 'math-clip a nil w))
X    ((and w (not (integerp w)))
X     (math-reject-arg w 'fixnump))
X    ((not (Math-num-integerp a))
X     (math-reject-arg a 'integerp))
X    ((< (or w (setq w calc-word-size)) 0)
X     (setq a (math-clip a (- w)))
X     (if (Math-natnum-lessp a (math-power-of-2 (- -1 w)))
X         a
X       (math-sub a (math-power-of-2 (- w)))))
X    ((Math-negp a)
X     (math-normalize (cons 'bigpos (math-binary-arg a w))))
X    ((and (integerp a) (< a 1000000))
X     (if (>= w 20)
X         a
X       (logand a (1- (lsh 1 w)))))
X    (t
X     (math-normalize
X      (cons 'bigpos
X        (math-clip-bignum (cdr (math-bignum-test (math-trunc a)))
X                  w)))))
)
(fset 'calcFunc-clip (symbol-function 'math-clip))
X
(defun math-clip-bignum (a w)   ; [l l]
X  (let ((q (math-div-bignum-digit a 512)))
X    (if (<= w 9)
X    (list (logand (cdr q)
X              (1- (lsh 1 w))))
X      (math-mul-bignum-digit (math-clip-bignum (math-norm-bignum (car q))
X                        (- w 9))
X                  512
X                  (cdr q))))
)
X
X
X
X
(defvar math-max-digits-cache nil)
(defun math-compute-max-digits (w r)
X  (let* ((pair (+ (* r 100000) w))
X     (res (assq pair math-max-digits-cache)))
X    (if res
X    (cdr res)
X      (let* ((calc-command-flags nil)
X         (digs (math-ceiling (math-div w (math-real-log2 r)))))
X    (setq math-max-digits-cache (cons (cons pair digs)
X                      math-max-digits-cache))
X    digs)))
)
X
(defvar math-log2-cache (list '(2 . 1)
X                  '(4 . 2)
X                  '(8 . 3)
X                  '(10 . (float 332193 -5))
X                  '(16 . 4)
X                  '(32 . 5)))
(defun math-real-log2 (x)   ;;; calc-internal-prec must be 6
X  (let ((res (assq x math-log2-cache)))
X    (if res
X    (cdr res)
X      (let* ((calc-symbolic-mode nil)
X         (calc-display-working-message nil)
X         (log (calcFunc-log x 2)))
X    (setq math-log2-cache (cons (cons x log) math-log2-cache))
X    log)))
)
X
(defconst math-radix-digits ["0" "1" "2" "3" "4" "5" "6" "7" "8" "9"
X                 "A" "B" "C" "D" "E" "F" "G" "H" "I" "J"
X                 "K" "L" "M" "N" "O" "P" "Q" "R" "S" "T"
X                 "U" "V" "W" "X" "Y" "Z"])
X
(defun math-format-radix (a)   ; [X S]
X  (if (< a calc-number-radix)
X      (if (< a 0)
X      (concat "-" (math-format-radix (- a)))
X    (math-format-radix-digit a))
X    (let ((s ""))
X      (while (> a 0)
X    (setq s (concat (math-format-radix-digit (% a calc-number-radix)) s)
X          a (/ a calc-number-radix)))
X      s))
)
X
(defconst math-binary-digits ["000" "001" "010" "011"
X                  "100" "101" "110" "111"])
(defun math-format-binary (a)   ; [X S]
X  (if (< a 8)
X      (if (< a 0)
X      (concat "-" (math-format-binary (- a)))
X    (math-format-radix a))
X    (let ((s ""))
X      (while (> a 7)
X    (setq s (concat (aref math-binary-digits (% a 8)) s)
X          a (/ a 8)))
X      (concat (math-format-radix a) s)))
)
X
(defun math-format-bignum-radix (a)   ; [X L]
X  (cond ((null a) "0")
X    ((and (null (cdr a))
X          (< (car a) calc-number-radix))
X     (math-format-radix-digit (car a)))
X    (t
X     (let ((q (math-div-bignum-digit a calc-number-radix)))
X       (concat (math-format-bignum-radix (math-norm-bignum (car q)))
X           (math-format-radix-digit (cdr q))))))
)
X
(defun math-format-bignum-binary (a)   ; [X L]
X  (cond ((null a) "0")
X    ((null (cdr a))
X     (math-format-binary (car a)))
X    (t
X     (let ((q (math-div-bignum-digit a 512)))
X       (concat (math-format-bignum-binary (math-norm-bignum (car q)))
X           (aref math-binary-digits (/ (cdr q) 64))
X           (aref math-binary-digits (% (/ (cdr q) 8) 8))
X           (aref math-binary-digits (% (cdr q) 8))))))
)
X
(defun math-format-bignum-octal (a)   ; [X L]
X  (cond ((null a) "0")
X    ((null (cdr a))
X     (math-format-radix (car a)))
X    (t
X     (let ((q (math-div-bignum-digit a 512)))
X       (concat (math-format-bignum-octal (math-norm-bignum (car q)))
X           (math-format-radix-digit (/ (cdr q) 64))
X           (math-format-radix-digit (% (/ (cdr q) 8) 8))
X           (math-format-radix-digit (% (cdr q) 8))))))
)
X
(defun math-format-bignum-hex (a)   ; [X L]
X  (cond ((null a) "0")
X    ((null (cdr a))
X     (math-format-radix (car a)))
X    (t
X     (let ((q (math-div-bignum-digit a 256)))
X       (concat (math-format-bignum-hex (math-norm-bignum (car q)))
X           (math-format-radix-digit (/ (cdr q) 16))
X           (math-format-radix-digit (% (cdr q) 16))))))
)
X
;;; Decompose into integer and fractional parts, without depending
;;; on calc-internal-prec.
(defun math-float-parts (a need-frac)    ; returns ( int frac fracdigs )
X  (if (>= (nth 2 a) 0)
X      (list (math-scale-rounding (nth 1 a) (nth 2 a)) '(float 0 0) 0)
X    (let* ((d (math-numdigs (nth 1 a)))
X       (n (- (nth 2 a))))
X      (if need-frac
X      (if (>= n d)
X          (list 0 a n)
X        (let ((qr (math-idivmod (nth 1 a) (math-scale-int 1 n))))
X          (list (car qr) (math-make-float (cdr qr) (- n)) n)))
X    (list (math-scale-rounding (nth 1 a) (nth 2 a))
X          '(float 0 0) 0))))
)
X
(defun math-format-radix-float (a prec)
X  (let ((fmt (car calc-float-format))
X    (figs (nth 1 calc-float-format))
X    (point calc-point-char)
X    (str nil))
X    (if (eq fmt 'fix)
X    (let* ((afigs (math-abs figs))
X           (fp (math-float-parts a (> afigs 0)))
X           (calc-internal-prec (+ 3 (max (nth 2 fp)
X                         (math-convert-radix-digits
X                          afigs t))))
X           (int (car fp))
X           (frac (math-round (math-mul (math-normalize (nth 1 fp))
X                       (math-radix-float-power afigs)))))
X      (if (not (and (math-zerop frac) (math-zerop int) (< figs 0)))
X          (let ((math-radix-explicit-format nil))
X        (let ((calc-group-digits nil))
X          (setq str (if (> afigs 0) (math-format-number frac) ""))
X          (if (< (length str) afigs)
X              (setq str (concat (make-string (- afigs (length str)) ?0)
X                    str))
X            (if (> (length str) afigs)
X            (setq str (substring str 1)
X                  int (math-add int 1))))
X          (setq str (concat (math-format-number int) point str)))
X        (if calc-group-digits
X            (setq str (math-group-float str))))
X        (setq figs 0))))
X    (or str
X    (let* ((prec calc-internal-prec)
X           (afigs (if (> figs 0)
X              figs
X            (max 1 (+ figs
X                  (1- (math-convert-radix-digits
X                       (max prec
X                        (math-numdigs (nth 1 a)))))))))
X           (calc-internal-prec (+ 3 (math-convert-radix-digits afigs t)))
X           (explo -1) (vlo (math-radix-float-power explo))
X           (exphi 1) (vhi (math-radix-float-power exphi))
X           expmid vmid eadj)
X      (setq a (math-normalize a))
X      (if (Math-zerop a)
X          (setq explo 0)
X        (if (math-lessp-float '(float 1 0) a)
X        (while (not (math-lessp-float a vhi))
X          (setq explo exphi vlo vhi
X            exphi (math-mul exphi 2)
X            vhi (math-radix-float-power exphi)))
X          (while (math-lessp-float a vlo)
X        (setq exphi explo vhi vlo
X              explo (math-mul explo 2)
X              vlo (math-radix-float-power explo))))
X        (while (not (eq (math-sub exphi explo) 1))
X          (setq expmid (math-div2 (math-add explo exphi))
X            vmid (math-radix-float-power expmid))
X          (if (math-lessp-float a vmid)
X          (setq exphi expmid vhi vmid)
X        (setq explo expmid vlo vmid)))
X        (setq a (math-div-float a vlo)))
X      (let* ((sc (math-round (math-mul a (math-radix-float-power
X                          (1- afigs)))))
X         (math-radix-explicit-format nil))
X        (let ((calc-group-digits nil))
X          (setq str (math-format-number sc))))
X      (if (> (length str) afigs)
X          (setq str (substring str 0 -1)
X            explo (1+ explo)))
X      (if (and (eq fmt 'float)
X           (math-lessp explo (+ (if (= figs 0)
X                        (1- (math-convert-radix-digits
X                         prec))
X                      afigs)
X                    calc-display-sci-high 1))
X           (math-lessp calc-display-sci-low explo))
X          (let ((dpos (1+ explo)))
X        (cond ((<= dpos 0)
X               (setq str (concat "0" point (make-string (- dpos) ?0)
X                     str)))
X              ((> dpos (length str))
X               (setq str (concat str (make-string (- dpos (length str))
X                              ?0) point)))
X              (t
X               (setq str (concat (substring str 0 dpos) point
X                     (substring str dpos)))))
X        (setq explo nil))
X        (setq eadj (if (eq fmt 'eng)
X               (min (math-mod explo 3) (length str))
X             0)
X          str (concat (substring str 0 (1+ eadj)) point
X                  (substring str (1+ eadj)))))
X      (setq pos (length str))
X      (while (eq (aref str (1- pos)) ?0) (setq pos (1- pos)))
X      (and explo (eq (aref str (1- pos)) ?.) (setq pos (1- pos)))
X      (setq str (substring str 0 pos))
X      (if calc-group-digits
X          (setq str (math-group-float str)))
X      (if explo
X          (let ((estr (let ((calc-number-radix 10)
X                (calc-group-digits nil))
X                (setq estr (math-format-number
X                    (math-sub explo eadj))))))
X        (setq str (if (or (memq calc-language '(math maple))
X                  (> calc-number-radix 14))
X                  (format "%s*%d.^%s" str calc-number-radix estr)
X                (format "%se%s" str estr)))))))
X    str)
)
X
(defun math-convert-radix-digits (n &optional to-dec)
X  (let ((key (cons n (cons to-dec calc-number-radix))))
X    (or (cdr (assoc key math-radix-digits-cache))
X    (let* ((calc-internal-prec 6)
X           (log (math-div (math-real-log2 calc-number-radix)
X                  '(float 332193 -5))))
X      (cdr (car (setq math-radix-digits-cache
X              (cons (cons key (math-ceiling (if to-dec
X                                (math-mul n log)
X                              (math-div n log))))
X                math-radix-digits-cache)))))))
)
(setq math-radix-digits-cache nil)
X
(defun math-radix-float-power (n)
X  (if (eq n 0)
X      '(float 1 0)
X    (or (and (eq calc-number-radix (car math-radix-float-cache-tag))
X         (<= calc-internal-prec (cdr math-radix-float-cache-tag)))
X    (setq math-radix-float-cache-tag (cons calc-number-radix
X                           calc-internal-prec)
X          math-radix-float-cache nil))
X    (math-normalize
X     (or (cdr (assoc n math-radix-float-cache))
X     (cdr (car (setq math-radix-float-cache
X             (cons (cons
X                n
X                (let ((calc-internal-prec
X                       (cdr math-radix-float-cache-tag)))
X                  (if (math-negp n)
X                      (math-div-float '(float 1 0)
X                              (math-radix-float-power
X                               (math-neg n)))
X                    (math-mul-float (math-sqr-float
X                             (math-radix-float-power
X                              (math-div2 n)))
X                            (if (math-evenp n)
X                            '(float 1 0)
X                              (math-float
X                               calc-number-radix))))))
X                   math-radix-float-cache)))))))
)
(setq math-radix-float-cache-tag nil)
X
SHAR_EOF
echo 'File calc-bin.el is complete' &&
chmod 0644 calc-bin.el ||
echo 'restore of calc-bin.el failed'
Wc_c="`wc -c < 'calc-bin.el'`"
test 24864 -eq "$Wc_c" ||
    echo 'calc-bin.el: original size 24864, current size' "$Wc_c"
rm -f _shar_wnt_.tmp
fi
# ============= calc-comb.el ==============
if test -f 'calc-comb.el' -a X"$1" != X"-c"; then
    echo 'x - skipping calc-comb.el (File already exists)'
    rm -f _shar_wnt_.tmp
else
> _shar_wnt_.tmp
echo 'x - extracting calc-comb.el (Text)'
sed 's/^X//' << 'SHAR_EOF' > 'calc-comb.el' &&
;; Calculator for GNU Emacs, part II [calc-comb.el]
;; Copyright (C) 1990, 1991 Free Software Foundation, Inc.
;; Written by Dave Gillespie, [e-mail].
X
;; This file is part of GNU Emacs.
X
;; GNU Emacs is distributed in the hope that it will be useful,
;; but WITHOUT ANY WARRANTY.  No author or distributor
;; accepts responsibility to anyone for the consequences of using it
;; or for whether it serves any particular purpose or works at all,
;; unless he says so in writing.  Refer to the GNU Emacs General Public
;; License for full details.
X
;; Everyone is granted permission to copy, modify and redistribute
;; GNU Emacs, but only under the conditions described in the
;; GNU Emacs General Public License.   A copy of this license is
;; supposed to have been given to you along with GNU Emacs so you
;; can know your rights and responsibilities.  It should be in a
;; file named COPYING.  Among other things, the copyright notice
;; and this notice must be preserved on all copies.
X
X
X
;; This file is autoloaded from calc-ext.el.
(require 'calc-ext)
X
(require 'calc-macs)
X
(defun calc-Need-calc-comb () nil)
X
X
;;; Combinatorics
X
(defun calc-gcd (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (calc-binary-op "gcd" 'calcFunc-gcd arg))
)
X
(defun calc-lcm (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (calc-binary-op "lcm" 'calcFunc-lcm arg))
)
X
(defun calc-extended-gcd ()
X  (interactive)
X  (calc-slow-wrapper
X   (calc-enter-result 2 "egcd" (cons 'calcFunc-egcd (calc-top-list-n 2))))
)
X
(defun calc-factorial (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (calc-unary-op "fact" 'calcFunc-fact arg))
)
X
(defun calc-gamma (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (calc-unary-op "gmma" 'calcFunc-gamma arg))
)
X
(defun calc-double-factorial (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (calc-unary-op "dfac" 'calcFunc-dfact arg))
)
X
(defun calc-choose (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (if (calc-is-hyperbolic)
X       (calc-binary-op "perm" 'calcFunc-perm arg)
X     (calc-binary-op "chos" 'calcFunc-choose arg)))
)
X
(defun calc-perm (arg)
X  (interactive "P")
X  (calc-hyperbolic-func)
X  (calc-choose arg)
)
X
(defvar calc-last-random-limit '(float 1 0))
(defun calc-random (n)
X  (interactive "P")
X  (calc-slow-wrapper
X   (if n
X       (calc-enter-result 0 "rand" (list 'calcFunc-random
X                     (calc-get-random-limit
X                      (prefix-numeric-value n))))
X     (calc-enter-result 1 "rand" (list 'calcFunc-random
X                       (calc-get-random-limit
X                    (calc-top-n 1))))))
)
X
(defun calc-get-random-limit (val)
X  (if (eq val 0)
X      calc-last-random-limit
X    (setq calc-last-random-limit val))
)
X
(defun calc-rrandom ()
X  (interactive)
X  (calc-slow-wrapper
X   (setq calc-last-random-limit '(float 1 0))
X   (calc-enter-result 0 "rand" (list 'calcFunc-random '(float 1 0))))
)
X
(defun calc-random-again (arg)
X  (interactive "p")
X  (calc-slow-wrapper
X   (while (>= (setq arg (1- arg)) 0)
X     (calc-enter-result 0 "rand" (list 'calcFunc-random
X                       calc-last-random-limit))))
)
X
(defun calc-shuffle (n)
X  (interactive "P")
X  (calc-slow-wrapper
X   (if n
X       (calc-enter-result 1 "shuf" (list 'calcFunc-shuffle
X                     (prefix-numeric-value n)
X                     (calc-get-random-limit
X                      (calc-top-n 1))))
X     (calc-enter-result 2 "shuf" (list 'calcFunc-shuffle
X                       (calc-top-n 1)
X                       (calc-get-random-limit
X                    (calc-top-n 2))))))
)
X
(defun calc-report-prime-test (res)
X  (cond ((eq (car res) t)
X     (calc-record-message "prim" "Prime (guaranteed)"))
X    ((eq (car res) nil)
X     (if (cdr res)
X         (if (eq (nth 1 res) 'unknown)
X         (calc-record-message
X          "prim" "Non-prime (factors unknown)")
X           (calc-record-message
X        "prim" "Non-prime (%s is a factor)"
X        (math-format-number (nth 1 res))))
X       (calc-record-message "prim" "Non-prime")))
X    (t
X     (calc-record-message
X      "prim" "Probably prime (%d iters; %s%% chance of error)"
X      (nth 1 res)
X      (let ((calc-float-format '(fix 2)))
X        (math-format-number (nth 2 res))))))
)
X
(defun calc-prime-test (iters)
X  (interactive "p")
X  (calc-slow-wrapper
X   (let* ((n (calc-top-n 1))
X      (res (math-prime-test n iters)))
X     (calc-report-prime-test res)))
)
X
(defun calc-next-prime (iters)
X  (interactive "p")
X  (calc-slow-wrapper
X   (let ((calc-verbose-nextprime t))
X     (if (calc-is-inverse)
X     (calc-enter-result 1 "prvp" (list 'calcFunc-prevprime
X                       (calc-top-n 1) (math-abs iters)))
X       (calc-enter-result 1 "nxtp" (list 'calcFunc-nextprime
X                     (calc-top-n 1) (math-abs iters))))))
)
X
(defun calc-prev-prime (iters)
X  (interactive "p")
X  (calc-invert-func)
X  (calc-next-prime iters)
)
X
(defun calc-prime-factors (iters)
X  (interactive "p")
X  (calc-slow-wrapper
X   (let ((res (calcFunc-prfac (calc-top-n 1))))
X     (if (not math-prime-factors-finished)
X     (calc-record-message "pfac" "Warning:  May not be fully factored"))
X     (calc-enter-result 1 "pfac" res)))
)
X
(defun calc-totient (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (calc-unary-op "phi" 'calcFunc-totient arg))
)
X
(defun calc-moebius (arg)
X  (interactive "P")
X  (calc-slow-wrapper
X   (calc-unary-op "mu" 'calcFunc-moebius arg))
)
X
X
X
X
X
(defun calcFunc-gcd (a b)
X  (if (Math-messy-integerp a)
X      (setq a (math-trunc a)))
X  (if (Math-messy-integerp b)
X      (setq b (math-trunc b)))
X  (cond ((and (Math-integerp a) (Math-integerp b))
X     (math-gcd a b))
X    ((Math-looks-negp a)
X     (calcFunc-gcd (math-neg a) b))
X    ((Math-looks-negp b)
X     (calcFunc-gcd a (math-neg b)))
X    ((Math-zerop a) b)
X    ((Math-zerop b) a)
X    ((and (Math-ratp a)
X          (Math-ratp b))
X     (math-make-frac (math-gcd (if (eq (car-safe a) 'frac) (nth 1 a) a)
X                   (if (eq (car-safe b) 'frac) (nth 1 b) b))
X             (calcFunc-lcm
X              (if (eq (car-safe a) 'frac) (nth 2 a) 1)
X              (if (eq (car-safe b) 'frac) (nth 2 b) 1))))
X    ((not (Math-integerp a))
X     (calc-record-why 'integerp a)
X     (list 'calcFunc-gcd a b))
X    (t
X     (calc-record-why 'integerp b)
X     (list 'calcFunc-gcd a b)))
)
X
(defun calcFunc-lcm (a b)
X  (let ((g (calcFunc-gcd a b)))
X    (if (Math-numberp g)
X    (math-div (math-mul a b) g)
X      (list 'calcFunc-lcm a b)))
)
X
(defun calcFunc-egcd (a b)   ; Knuth section 4.5.2
X  (cond
X   ((not (Math-integerp a))
X    (if (Math-messy-integerp a)
X    (calcFunc-egcd (math-trunc a) b)
X      (calc-record-why 'integerp a)
X      (list 'calcFunc-egcd a b)))
X   ((not (Math-integerp b))
X    (if (Math-messy-integerp b)
X    (calcFunc-egcd a (math-trunc b))
X      (calc-record-why 'integerp b)
X      (list 'calcFunc-egcd a b)))
X   (t
X    (let ((u1 1) (u2 0) (u3 a)
X      (v1 0) (v2 1) (v3 b)
X      t1 t2 q)
X      (while (not (eq v3 0))
X    (setq q (math-idivmod u3 v3)
X          t1 (math-sub u1 (math-mul v1 (car q)))
X          t2 (math-sub u2 (math-mul v2 (car q)))
X          u1 v1  u2 v2  u3 v3
X          v1 t1  v2 t2  v3 (cdr q)))
X      (list 'vec u3 u1 u2))))
)
X
X
;;; Factorial and related functions.
X
(defun calcFunc-fact (n)   ; [I I] [F F] [Public]
X  (let (temp)
X    (cond ((Math-integer-negp n)
X       (if calc-infinite-mode
X           '(var uinf var-uinf)
X         (math-reject-arg n 'range)))
X      ((integerp n)
X       (if (<= n 20)
X           (aref '[1 1 2 6 24 120 720 5040 40320 362880
X             (bigpos 800 628 3) (bigpos 800 916 39)
X             (bigpos 600 1 479) (bigpos 800 20 227 6)
X             (bigpos 200 291 178 87) (bigpos 0 368 674 307 1)
X             (bigpos 0 888 789 922 20) (bigpos 0 96 428 687 355)
X             (bigpos 0 728 705 373 402 6)
X             (bigpos 0 832 408 100 645 121)
X             (bigpos 0 640 176 8 902 432 2)] n)
X         (math-factorial-iter (1- n) 2 1)))
X      ((and (math-messy-integerp n)
X        (Math-lessp n 100))
X       (math-inexact-result)
X       (setq temp (math-trunc n))
X       (if (>= temp 0)
X           (if (<= temp 20)
X           (math-float (calcFunc-fact temp))
X         (math-with-extra-prec 1
X           (math-factorial-iter (1- temp) 2 '(float 1 0))))
X         (math-reject-arg n 'range)))
X      ((math-numberp n)
X       (let* ((q (math-quarter-integer n))
X          (tn (and q (Math-lessp n 1000) (Math-lessp -1000 n)
X               (1+ (math-floor n)))))
X         (cond ((and tn (= q 2)
X             (or calc-symbolic-mode (< (math-abs tn) 20)))
X            (let ((q (if (< tn 0)
X                 (math-div
X                  (math-pow -2 (- tn))
X                  (math-double-factorial-iter (* -2 tn) 3 1 2))
X                   (math-div
X                (math-double-factorial-iter (* 2 tn) 3 1 2)
X                (math-pow 2 tn)))))
X              (math-mul q (if calc-symbolic-mode
X                      (list 'calcFunc-sqrt '(var pi var-pi))
X                    (math-sqrt-pi)))))
X           ((and tn (>= tn 0) (< tn 20)
X             (memq q '(1 3)))
X            (math-inexact-result)
X            (math-div
X             (math-mul (math-double-factorial-iter (* 4 tn) q 1 4)
X                   (if (= q 1) (math-gamma-1q) (math-gamma-3q)))
X             (math-pow 4 tn)))
X           (t
X            (math-inexact-result)
X            (math-with-extra-prec 3
X              (math-gammap1-raw (math-float n)))))))
X      ((equal n '(var inf var-inf)) n)
X      (t (calc-record-why 'numberp n)
X         (list 'calcFunc-fact n))))
)
X
(math-defcache math-gamma-1q nil
X  (math-with-extra-prec 3
X    (math-gammap1-raw '(float -75 -2))))
X
(math-defcache math-gamma-3q nil
X  (math-with-extra-prec 3
X    (math-gammap1-raw '(float -25 -2))))
X
(defun math-factorial-iter (count n f)
X  (if (= (% n 5) 1)
X      (math-working (format "factorial(%d)" (1- n)) f))
X  (if (> count 0)
X      (math-factorial-iter (1- count) (1+ n) (math-mul n f))
X    f)
)
X
(defun calcFunc-dfact (n)   ; [I I] [F F] [Public]
X  (cond ((Math-integer-negp n)
X     (if (math-oddp n)
X         (if (eq n -1)
X         1
X           (math-div (if (eq (math-mod n 4) 3) 1 -1)
X             (calcFunc-dfact (math-sub -2 n))))
X       (list 'calcFunc-dfact n)))
X    ((Math-zerop n) 1)
X    ((integerp n) (math-double-factorial-iter n (+ 2 (% n 2)) 1 2))
X    ((math-messy-integerp n)
X     (let ((temp (math-trunc n)))
X       (math-inexact-result)
X       (if (natnump temp)
X           (if (Math-lessp temp 200)
X           (math-with-extra-prec 1
X             (math-double-factorial-iter temp (+ 2 (% temp 2))
X                         '(float 1 0) 2))
X         (let* ((half (math-div2 temp))
X            (even (math-mul (math-pow 2 half)
X                    (calcFunc-fact (math-float half)))))
X           (if (math-evenp temp)
X               even
X             (math-div (calcFunc-fact n) even))))
X         (list 'calcFunc-dfact max))))
X    ((equal n '(var inf var-inf)) n)
X    (t (calc-record-why 'natnump n)
X       (list 'calcFunc-dfact n)))
)
X
(defun math-double-factorial-iter (max n f step)
X  (if (< (% n 12) step)
X      (math-working (format "dfact(%d)" (- n step)) f))
X  (if (<= n max)
X      (math-double-factorial-iter max (+ n step) (math-mul n f) step)
X    f)
)
X
(defun calcFunc-perm (n m)   ; [I I I] [F F F] [Public]
X  (cond ((and (integerp n) (integerp m) (<= m n) (>= m 0))
X     (math-factorial-iter m (1+ (- n m)) 1))
X    ((or (not (math-num-integerp n))
X         (and (math-messy-integerp n) (Math-lessp 100 n))
X         (not (math-num-integerp m))
X         (and (math-messy-integerp m) (Math-lessp 100 m)))
X     (or (math-realp n) (equal n '(var inf var-inf))
X         (math-reject-arg n 'realp))
X     (or (math-realp m) (equal m '(var inf var-inf))
X         (math-reject-arg m 'realp))
X     (and (math-num-integerp n) (math-negp n) (math-reject-arg n 'range))
X     (and (math-num-integerp m) (math-negp m) (math-reject-arg m 'range))
X     (math-div (calcFunc-fact n) (calcFunc-fact (math-sub n m))))
X    (t
X     (let ((tn (math-trunc n))
X           (tm (math-trunc m)))
X       (math-inexact-result)
X       (or (integerp tn) (math-reject-arg tn 'fixnump))
X       (or (integerp tm) (math-reject-arg tm 'fixnump))
X       (or (and (<= tm tn) (>= tm 0)) (math-reject-arg tm 'range))
X       (math-with-extra-prec 1
X         (math-factorial-iter tm (1+ (- tn tm)) '(float 1 0))))))
)
X
(defun calcFunc-choose (n m)   ; [I I I] [F F F] [Public]
X  (cond ((and (integerp n) (integerp m) (<= m n) (>= m 0))
X     (if (> m (/ n 2))
X         (math-choose-iter (- n m) n 1 1)
X       (math-choose-iter m n 1 1)))
X    ((not (math-realp n))
X     (math-reject-arg n 'realp))
X    ((not (math-realp m))
X     (math-reject-arg m 'realp))
X    ((not (math-num-integerp m))
X     (if (and (math-num-integerp n) (math-negp n))
X         (list 'calcFunc-choose n m)
X       (math-div (calcFunc-fact (math-float n))
X             (math-mul (calcFunc-fact m)
X                   (calcFunc-fact (math-sub n m))))))
X    ((math-negp m) 0)
X    ((math-negp n)
X     (let ((val (calcFunc-choose (math-add (math-add n m) -1) m)))
X       (if (math-evenp (math-trunc m))
X           val
X         (math-neg val))))
X    ((and (math-num-integerp n)
X          (Math-lessp n m))
X     0)
X    (t
X     (math-inexact-result)
X     (let ((tm (math-trunc m)))
X       (or (integerp tm) (math-reject-arg tm 'fixnump))
X       (if (> tm 100)
X           (math-div (calcFunc-fact (math-float n))
X             (math-mul (calcFunc-fact (math-float m))
X                   (calcFunc-fact (math-float
X                           (math-sub n m)))))
X         (math-with-extra-prec 1
X           (math-choose-float-iter tm n 1 1))))))
)
X
(defun math-choose-iter (m n i c)
X  (if (and (= (% i 5) 1) (> i 5))
X      (math-working (format "choose(%d)" (1- i)) c))
X  (if (<= i m)
X      (math-choose-iter m (1- n) (1+ i)
X            (math-quotient (math-mul c n) i))
X    c)
)
X
(defun math-choose-float-iter (count n i c)
X  (if (= (% i 5) 1)
X      (math-working (format "choose(%d)" (1- i)) c))
X  (if (> count 0)
X      (math-choose-float-iter (1- count) (math-sub n 1) (1+ i)
X                  (math-div (math-mul c n) i))
X    c)
)
X
X
;;; Stirling numbers.
X
(defun calcFunc-stir1 (n m)
X  (math-stirling-number n m 1)
)
X
(defun calcFunc-stir2 (n m)
X  (math-stirling-number n m 0)
)
X
(defun math-stirling-number (n m k)
X  (or (math-num-natnump n) (math-reject-arg n 'natnump))
X  (or (math-num-natnump m) (math-reject-arg m 'natnump))
X  (if (consp n) (setq n (math-trunc n)))
X  (or (integerp n) (math-reject-arg n 'fixnump))
X  (if (consp m) (setq m (math-trunc m)))
X  (or (integerp m) (math-reject-arg m 'fixnump))
X  (if (< n m)
X      0
X    (let ((cache (aref math-stirling-cache k)))
X      (while (<= (length cache) n)
X    (let ((i (1- (length cache)))
X          row)
X      (setq cache (vconcat cache (make-vector (length cache) nil)))
X      (aset math-stirling-cache k cache)
X      (while (< (setq i (1+ i)) (length cache))
X        (aset cache i (setq row (make-vector (1+ i) nil)))
X        (aset row 0 0)
X        (aset row i 1))))
X      (if (= k 1)
X      (math-stirling-1 n m)
X    (math-stirling-2 n m))))
)
(setq math-stirling-cache (vector [[1]] [[1]]))
X
(defun math-stirling-1 (n m)
X  (or (aref (aref cache n) m)
X      (aset (aref cache n) m
X        (math-add (math-stirling-1 (1- n) (1- m))
X              (math-mul (- 1 n) (math-stirling-1 (1- n) m)))))
)
X
(defun math-stirling-2 (n m)
X  (or (aref (aref cache n) m)
X      (aset (aref cache n) m
X        (math-add (math-stirling-2 (1- n) (1- m))
X              (math-mul m (math-stirling-2 (1- n) m)))))
)
X
X
;;; Produce a random 10-bit integer, with (random) if no seed provided,
;;; or else with Numerical Recipes algorithm ran3 / Knuth 3.2.2-A.
(defun math-init-random-base ()
X  (if (and (boundp 'var-RandSeed) var-RandSeed)
X      (if (eq (car-safe var-RandSeed) 'vec)
X      nil
X    (if (Math-integerp var-RandSeed)
X        (let* ((seed (math-sub 161803 var-RandSeed))
X           (mj (1+ (math-mod seed '(bigpos 0 0 1))))
X           (mk (1+ (math-mod (math-quotient seed '(bigpos 0 0 1))
X                     '(bigpos 0 0 1))))
X           (i 0))
X          (setq math-random-table (cons 'vec (make-list 55 mj)))
X          (while (<= (setq i (1+ i)) 54)
X        (let* ((ii (% (* i 21) 55))
X               (p (nthcdr ii math-random-table)))
X          (setcar p mk)
X          (setq mk (- mj mk)
X            mj (car p)))))
X      (math-reject-arg var-RandSeed "*RandSeed must be an integer"))
X    (setq var-RandSeed (list 'vec var-RandSeed)
X          math-random-ptr1 math-random-table
X          math-random-cache nil
X          math-random-ptr2 (nthcdr 31 math-random-table))
X    (let ((i 200))
X      (while (> (setq i (1- i)) 0)
X        (math-random-base))))
X    (random t)
X    (setq var-RandSeed nil
X      math-random-cache nil
X      i 0
X      math-random-shift -4)  ; assume RAND_MAX >= 16383
X    ;; This exercises the random number generator and also helps
X    ;; deduce a better value for RAND_MAX.
X    (while (< (setq i (1+ i)) 30)
X      (if (> (lsh (math-abs (random)) math-random-shift) 4095)
X      (setq math-random-shift (1- math-random-shift)))))
X  (setq math-last-RandSeed var-RandSeed
X    math-gaussian-cache nil)
)
X
(defun math-random-base ()
X  (if var-RandSeed
X      (progn
X    (setq math-random-ptr1 (or (cdr math-random-ptr1)
X                   (cdr math-random-table))
X          math-random-ptr2 (or (cdr math-random-ptr2)
X                   (cdr math-random-table)))
X    (logand (lsh (setcar math-random-ptr1
X                 (logand (- (car math-random-ptr1)
X                    (car math-random-ptr2)) 524287))
X             -6) 1023))
X    (logand (lsh (random) math-random-shift) 1023))
)
(setq math-random-table nil)
(setq math-last-RandSeed nil)
(setq math-random-ptr1 nil)
(setq math-random-ptr2 nil)
(setq math-random-shift nil)
X
X
;;; Produce a random digit in the range 0..999.
;;; Avoid various pitfalls that may lurk in the built-in (random) function!
;;; Shuffling algorithm from Numerical Recipes, section 7.1.
(defun math-random-digit ()
X  (let (i)
X    (or (and (boundp 'var-RandSeed) (eq var-RandSeed math-last-RandSeed))
X    (math-init-random-base))
X    (or math-random-cache
X    (progn
X      (setq math-random-last (math-random-base)
X        math-random-cache (make-vector 13 nil)
X        i -1)
X      (while (< (setq i (1+ i)) 13)
X        (aset math-random-cache i (math-random-base)))))
X    (while (progn
X         (setq i (/ math-random-last 79)   ; 0 <= i < 13
X           math-random-last (aref math-random-cache i))
X         (aset math-random-cache i (math-random-base))
X         (>= math-random-last 1000)))
X    math-random-last)
)
(setq math-random-cache nil)
X
;;; Produce an N-digit random integer.
(defun math-random-digits (n)
X  (cond ((<= n 6)
X     (math-scale-right (+ (* (math-random-digit) 1000) (math-random-digit))
X               (- 6 n)))
X    (t (let* ((slop (% (- 900003 n) 3))
X          (i (/ (+ n slop) 3))
X          (digs nil))
X         (while (> i 0)
X           (setq digs (cons (math-random-digit) digs)
X             i (1- i)))
X         (math-normalize (math-scale-right (cons 'bigpos digs)
X                           slop)))))
)
X
;;; Produce a uniformly-distributed random float 0 <= N < 1.
(defun math-random-float ()
X  (math-make-float (math-random-digits calc-internal-prec)
X           (- calc-internal-prec))
)
X
;;; Produce a Gaussian-distributed random float with mean=0, sigma=1.
(defun math-gaussian-float ()
X  (math-with-extra-prec 2
X    (if (and math-gaussian-cache
X         (= (car math-gaussian-cache) calc-internal-prec))
X    (prog1
X        (cdr math-gaussian-cache)
X      (setq math-gaussian-cache nil))
X      (let* ((v1 (math-add (math-mul (math-random-float) 2) -1))
X         (v2 (math-add (math-mul (math-random-float) 2) -1))
X         (r (math-add (math-sqr v1) (math-sqr v2))))
X    (while (or (not (Math-lessp r 1)) (math-zerop r))
X      (setq v1 (math-add (math-mul (math-random-float) 2) -1)
X        v2 (math-add (math-mul (math-random-float) 2) -1)
X        r (math-add (math-sqr v1) (math-sqr v2))))
X    (let ((fac (math-sqrt (math-mul (math-div (calcFunc-ln r) r) -2))))
X      (setq math-gaussian-cache (cons calc-internal-prec
X                      (math-mul v1 fac)))
X      (math-mul v2 fac)))))
)
(setq math-gaussian-cache nil)
X
;;; Produce a random integer or real 0 <= N < MAX.
(defun calcFunc-random (max)
X  (cond ((Math-zerop max)
X     (math-gaussian-float))
X    ((Math-integerp max)
X     (let* ((digs (math-numdigs max))
X        (r (math-random-digits (+ digs 3))))
X       (math-mod r max)))
X    ((Math-realp max)
X     (math-mul (math-random-float) max))
X    ((and (eq (car max) 'intv) (math-constp max)
X          (Math-lessp (nth 2 max) (nth 3 max)))
X     (if (math-floatp max)
X         (let ((val (math-add (math-mul (math-random-float)
X                        (math-sub (nth 3 max) (nth 2 max)))
X                  (nth 2 max))))
X           (if (or (and (memq (nth 1 max) '(0 1))      ; almost not worth
X                (Math-equal val (nth 2 max)))  ;   checking!
X               (and (memq (nth 1 max) '(0 2))
X                (Math-equal val (nth 3 max))))
X           (calcFunc-random max)
X         val))
X       (let ((lo (if (memq (nth 1 max) '(0 1))
X             (math-add (nth 2 max) 1) (nth 2 max)))
X         (hi (if (memq (nth 1 max) '(1 3))
X             (math-add (nth 3 max) 1) (nth 3 max))))
X         (if (Math-lessp lo hi)
X         (math-add (calcFunc-random (math-sub hi lo)) lo)
X           (math-reject-arg max "*Empty interval")))))
X    ((eq (car max) 'vec)
X     (if (cdr max)
X         (nth (1+ (calcFunc-random (1- (length max)))) max)
X       (math-reject-arg max "*Empty list")))
X    ((and (eq (car max) 'sdev) (math-constp max) (Math-realp (nth 1 max)))
X     (math-add (math-mul (math-gaussian-float) (nth 2 max)) (nth 1 max)))
X    (t (math-reject-arg max 'realp)))
)
X
;;; Choose N objects at random from the set MAX without duplicates.
(defun calcFunc-shuffle (n &optional max)
X  (or max (setq max n n -1))
X  (or (and (Math-num-integerp n)
X       (or (natnump (setq n (math-trunc n))) (eq n -1)))
X      (math-reject-arg n 'integerp))
X  (cond ((or (math-zerop max)
X         (math-floatp max)
X         (eq (car-safe max) 'sdev))
X     (if (< n 0)
X         (math-reject-arg n 'natnump)
X       (math-simple-shuffle n max)))
X    ((and (<= n 1) (>= n 0))
X     (math-simple-shuffle n max))
X    ((and (eq (car-safe max) 'intv) (math-constp max))
X     (let ((num (math-add (math-sub (nth 3 max) (nth 2 max))
X                  (cdr (assq (nth 1 max)
X                     '((0 . -1) (1 . 0)
X                       (2 . 0) (3 . 1))))))
X           (min (math-add (nth 2 max) (if (memq (nth 1 max) '(0 1))
X                          1 0))))
X       (if (< n 0) (setq n num))
X       (or (math-posp num) (math-reject-arg max 'range))
X       (and (Math-lessp num n) (math-reject-arg n 'range))
X       (if (Math-lessp n (math-quotient num 3))
X           (math-simple-shuffle n max)
X         (if (> (* n 4) (* num 3))
X         (math-add (math-sub min 1)
X               (math-shuffle-list n num (calcFunc-index num)))
X           (let ((tot 0)
X             (m 0)
X             (vec nil))
X         (while (< m n)
X           (if (< (calcFunc-random (- num tot)) (- n m))
X               (setq vec (cons (math-add min tot) vec)
X                 m (1+ m)))
X           (setq tot (1+ tot)))
X         (math-shuffle-list n n (cons 'vec vec)))))))
X    ((eq (car-safe max) 'vec)
X     (let ((size (1- (length max))))
X       (if (< n 0) (setq n size))
X       (if (and (> n (/ size 2)) (<= n size))
X           (math-shuffle-list n size (copy-sequence max))
X         (let* ((vals (calcFunc-shuffle
X               n (list 'intv 3 1 (1- (length max)))))
X            (p vals))
X           (while (setq p (cdr p))
X         (setcar p (nth (car p) max)))
X           vals))))
X    ((math-integerp max)
X     (if (math-posp max)
X         (calcFunc-shuffle n (list 'intv 2 0 max))
X       (calcFunc-shuffle n (list 'intv 1 max 0))))
X    (t (math-reject-arg max 'realp)))
)
X
(defun math-simple-shuffle (n max)
X  (let ((vec nil)
X    val)
X    (while (>= (setq n (1- n)) 0)
X      (while (math-member (setq val (calcFunc-random max)) vec))
X      (setq vec (cons val vec)))
X    (cons 'vec vec))
)
X
(defun math-shuffle-list (n size vec)
X  (let ((j size)
X    k temp
X    (p vec))
X    (while (cdr (setq p (cdr p)))
X      (setq k (calcFunc-random j)
X        j (1- j)
X        temp (nth k p))
X      (setcar (nthcdr k p) (car p))
X      (setcar p temp))
X    (cons 'vec (nthcdr (- size n -1) vec)))
)
X
(defun math-member (x list)
X  (while (and list (not (equal x (car list))))
X    (setq list (cdr list)))
X  list
)
X
X
;;; Check if the integer N is prime.  [X I]
;;; Return (nil) if non-prime,
;;;        (nil N) if non-prime with known factor N,
;;;        (nil unknown) if non-prime with no known factors,
;;;        (t) if prime,
;;;        (maybe N P) if probably prime (after N iters with probability P%)
(defun math-prime-test (n iters)
X  (if (and (Math-vectorp n) (cdr n))
X      (setq n (nth (1- (length n)) n)))
X  (if (Math-messy-integerp n)
X      (setq n (math-trunc n)))
X  (let ((res))
X    (while (> iters 0)
X      (setq res
X        (cond ((and (integerp n) (<= n 5003))
X           (list (= (math-next-small-prime n) n)))
X          ((not (Math-integerp n))
X           (error "Argument must be an integer"))
X          ((Math-integer-negp n)
X           '(nil))
X          ((Math-natnum-lessp n '(bigpos 0 0 8))
X           (setq n (math-fixnum n))
X           (let ((i -1) v)
X             (while (and (> (% n (setq v (aref math-primes-table
X                               (setq i (1+ i)))))
X                    0)
X                 (< (* v v) n)))
X             (if (= (% n v) 0)
X             (list nil v)
X               '(t))))
X          ((not (equal n (car math-prime-test-cache)))
X           (cond ((= (% (nth 1 n) 2) 0) '(nil 2))
X             ((= (% (nth 1 n) 5) 0) '(nil 5))
X             (t (let ((dig (cdr n)) (sum 0))
X                  (while dig
X                (if (cdr dig)
X                    (setq sum (% (+ (+ sum (car dig))
X                            (* (nth 1 dig) 1000))
X                         111111)
X                      dig (cdr (cdr dig)))
X                  (setq sum (% (+ sum (car dig)) 111111)
X                    dig nil)))
X                  (cond ((= (% sum 3) 0) '(nil 3))
X                    ((= (% sum 7) 0) '(nil 7))
X                    ((= (% sum 11) 0) '(nil 11))
X                    ((= (% sum 13) 0) '(nil 13))
X                    ((= (% sum 37) 0) '(nil 37))
X                    (t
X                     (setq math-prime-test-cache-k 1
X                       math-prime-test-cache-q
X                       (math-div2 n)
X                       math-prime-test-cache-nm1
X                       (math-add n -1))
X                     (while (math-evenp
X                         math-prime-test-cache-q)
X                       (setq math-prime-test-cache-k
X                         (1+ math-prime-test-cache-k)
X                         math-prime-test-cache-q
X                         (math-div2
X                          math-prime-test-cache-q)))
X                     (setq iters (1+ iters))
X                     (list 'maybe
X                       0
X                       (math-sub
X                        100
X                        (math-div
X                         '(float 232 0)
X                         (math-numdigs n))))))))))
X          ((not (eq (car (nth 1 math-prime-test-cache)) 'maybe))
X           (nth 1 math-prime-test-cache))
X          (t   ; Fermat step
X           (let* ((x (math-add (calcFunc-random (math-add n -2)) 2))
X              (y (math-pow-mod x math-prime-test-cache-q n))
X              (j 0))
X             (while (and (not (eq y 1))
X                 (not (equal y math-prime-test-cache-nm1))
X                 (< (setq j (1+ j)) math-prime-test-cache-k))
X               (setq y (math-mod (math-mul y y) n)))
X             (if (or (equal y math-prime-test-cache-nm1)
X                 (and (eq y 1) (eq j 0)))
X             (list 'maybe
X                   (1+ (nth 1 (nth 1 math-prime-test-cache)))
X                   (math-mul (nth 2 (nth 1 math-prime-test-cache))
X                     '(float 25 -2)))
X               '(nil unknown))))))
X      (setq math-prime-test-cache (list n res)
X        iters (if (eq (car res) 'maybe)
X              (1- iters)
X            0)))
X    res)
)
(defvar math-prime-test-cache '(-1))
X
(defun calcFunc-prime (n &optional iters)
X  (or (math-num-integerp n) (math-reject-arg n 'integerp))
X  (or (not iters) (math-num-integerp iters) (math-reject-arg iters 'integerp))
X  (if (car (math-prime-test (math-trunc n) (math-trunc (or iters 1))))
X      1
X    0)
)
X
;;; Theory: summing base-10^6 digits modulo 111111 is "casting out 999999s".
;;; Initial probability that N is prime is 1/ln(N) = log10(e)/log10(N).
;;; After culling [2,3,5,7,11,13,37], probability of primality is 5.36 x more.
;;; Initial reported probability of non-primality is thus 100% - this.
;;; Each Fermat step multiplies this probability by 25%.
;;; The Fermat step is algorithm P from Knuth section 4.5.4.
X
X
(defun calcFunc-prfac (n)
X  (setq math-prime-factors-finished t)
X  (if (Math-messy-integerp n)
X      (setq n (math-trunc n)))
X  (if (Math-natnump n)
X      (if (Math-natnum-lessp 2 n)
X      (let (factors res p (i 0))
X        (while (and (not (eq n 1))
X            (< i (length math-primes-table)))
X          (setq p (aref math-primes-table i))
X          (while (eq (cdr (setq res (cond ((eq n p) (cons 1 0))
X                          ((eq n 1) (cons 0 1))
X                          ((consp n) (math-idivmod n p))
X                          (t (cons (/ n p) (% n p))))))
X             0)
X        (math-working "factor" p)
X        (setq factors (nconc factors (list p))
X              n (car res)))
X          (or (eq n 1)
X          (Math-natnum-lessp p (car res))
X          (setq factors (nconc factors (list n))
X            n 1))
X          (setq i (1+ i)))
X        (or (setq math-prime-factors-finished (eq n 1))
X        (setq factors (nconc factors (list n))))
X        (cons 'vec factors))
X    (list 'vec n))
X    (if (Math-integerp n)
X    (if (eq n -1)
X        (list 'vec n)
X      (cons 'vec (cons -1 (cdr (calcFunc-prfac (math-neg n))))))
X      (calc-record-why 'integerp n)
X      (list 'calcFunc-prfac n)))
)
X
(defun calcFunc-totient (n)
X  (if (Math-messy-integerp n)
X      (setq n (math-trunc n)))
X  (if (Math-natnump n)
X      (if (Math-natnum-lessp n 2)
X      (if (Math-negp n)
X          (calcFunc-totient (math-abs n))
X        n)
X    (let ((factors (cdr (calcFunc-prfac n)))
X          p)
X      (if math-prime-factors-finished
X          (progn
X        (while factors
X          (setq p (car factors)
X            n (math-mul (math-div n p) (math-add p -1)))
X          (while (equal p (car factors))
X            (setq factors (cdr factors))))
X        n)
X        (calc-record-why "*Number too big to factor" n)
X        (list 'calcFunc-totient n))))
X    (calc-record-why 'natnump n)
X    (list 'calcFunc-totient n))
)
X
(defun calcFunc-moebius (n)
X  (if (Math-messy-integerp n)
X      (setq n (math-trunc n)))
X  (if (and (Math-natnump n) (not (eq n 0)))
X      (if (Math-natnum-lessp n 2)
X      (if (Math-negp n)
X          (calcFunc-moebius (math-abs n))
X        1)
X    (let ((factors (cdr (calcFunc-prfac n)))
X          (mu 1))
X      (if math-prime-factors-finished
X          (progn
X        (while factors
X          (setq mu (if (equal (car factors) (nth 1 factors))
X                   0 (math-neg mu))
X            factors (cdr factors)))
X        mu)
X        (calc-record-why "Number too big to factor" n)
X        (list 'calcFunc-moebius n))))
X    (calc-record-why 'posintp n)
X    (list 'calcFunc-moebius n))
)
X
X
(defun calcFunc-nextprime (n &optional iters)
X  (if (Math-integerp n)
X      (if (Math-integer-negp n)
X      2
X    (if (and (integerp n) (< n 5003))
X        (math-next-small-prime (1+ n))
X      (if (math-evenp n)
X          (setq n (math-add n -1)))
X      (let (res)
X        (while (not (car (setq res (math-prime-test
X                    (setq n (math-add n 2))
X                    (or iters 1))))))
X        (if (and calc-verbose-nextprime
X             (eq (car res) 'maybe))
X        (calc-report-prime-test res)))
X      n))
X    (if (Math-realp n)
X    (calcFunc-nextprime (math-trunc n) iters)
X      (math-reject-arg n 'integerp)))
)
(setq calc-verbose-nextprime nil)
X
(defun calcFunc-prevprime (n &optional iters)
X  (if (Math-integerp n)
X      (if (Math-lessp n 4)
X      2
X    (if (math-evenp n)
X        (setq n (math-add n 1)))
X    (let (res)
X      (while (not (car (setq res (math-prime-test
X                      (setq n (math-add n -2))
X                      (or iters 1))))))
X      (if (and calc-verbose-nextprime
X           (eq (car res) 'maybe))
X          (calc-report-prime-test res)))
X    n)
X    (if (Math-realp n)
X    (calcFunc-prevprime (math-ceiling n) iters)
X      (math-reject-arg n 'integerp)))
)
X
(defun math-next-small-prime (n)
X  (if (and (integerp n) (> n 2))
X      (let ((lo -1)
X        (hi (length math-primes-table))
X        mid)
X    (while (> (- hi lo) 1)
X      (if (> n (aref math-primes-table
X             (setq mid (ash (+ lo hi) -1))))
X          (setq lo mid)
X        (setq hi mid)))
X    (aref math-primes-table hi))
X    2)
)
X
(defconst math-primes-table
X  [2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89
X     97 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181
X     191 193 197 199 211 223 227 229 233 239 241 251 257 263 269 271 277
X     281 283 293 307 311 313 317 331 337 347 349 353 359 367 373 379 383
X     389 397 401 409 419 421 431 433 439 443 449 457 461 463 467 479 487
X     491 499 503 509 521 523 541 547 557 563 569 571 577 587 593 599 601
X     607 613 617 619 631 641 643 647 653 659 661 673 677 683 691 701 709
X     719 727 733 739 743 751 757 761 769 773 787 797 809 811 821 823 827
X     829 839 853 857 859 863 877 881 883 887 907 911 919 929 937 941 947
X     953 967 971 977 983 991 997 1009 1013 1019 1021 1031 1033 1039 1049
X     1051 1061 1063 1069 1087 1091 1093 1097 1103 1109 1117 1123 1129 1151
X     1153 1163 1171 1181 1187 1193 1201 1213 1217 1223 1229 1231 1237 1249
X     1259 1277 1279 1283 1289 1291 1297 1301 1303 1307 1319 1321 1327 1361
X     1367 1373 1381 1399 1409 1423 1427 1429 1433 1439 1447 1451 1453 1459
X     1471 1481 1483 1487 1489 1493 1499 1511 1523 1531 1543 1549 1553 1559
X     1567 1571 1579 1583 1597 1601 1607 1609 1613 1619 1621 1627 1637 1657
X     1663 1667 1669 1693 1697 1699 1709 1721 1723 1733 1741 1747 1753 1759
X     1777 1783 1787 1789 1801 1811 1823 1831 1847 1861 1867 1871 1873 1877
X     1879 1889 1901 1907 1913 1931 1933 1949 1951 1973 1979 1987 1993 1997
X     1999 2003 2011 2017 2027 2029 2039 2053 2063 2069 2081 2083 2087 2089
X     2099 2111 2113 2129 2131 2137 2141 2143 2153 2161 2179 2203 2207 2213
X     2221 2237 2239 2243 2251 2267 2269 2273 2281 2287 2293 2297 2309 2311
X     2333 2339 2341 2347 2351 2357 2371 2377 2381 2383 2389 2393 2399 2411
X     2417 2423 2437 2441 2447 2459 2467 2473 2477 2503 2521 2531 2539 2543
X     2549 2551 2557 2579 2591 2593 2609 2617 2621 2633 2647 2657 2659 2663
X     2671 2677 2683 2687 2689 2693 2699 2707 2711 2713 2719 2729 2731 2741
X     2749 2753 2767 2777 2789 2791 2797 2801 2803 2819 2833 2837 2843 2851
X     2857 2861 2879 2887 2897 2903 2909 2917 2927 2939 2953 2957 2963 2969
X     2971 2999 3001 3011 3019 3023 3037 3041 3049 3061 3067 3079 3083 3089
X     3109 3119 3121 3137 3163 3167 3169 3181 3187 3191 3203 3209 3217 3221
X     3229 3251 3253 3257 3259 3271 3299 3301 3307 3313 3319 3323 3329 3331
X     3343 3347 3359 3361 3371 3373 3389 3391 3407 3413 3433 3449 3457 3461
X     3463 3467 3469 3491 3499 3511 3517 3527 3529 3533 3539 3541 3547 3557
X     3559 3571 3581 3583 3593 3607 3613 3617 3623 3631 3637 3643 3659 3671
X     3673 3677 3691 3697 3701 3709 3719 3727 3733 3739 3761 3767 3769 3779
SHAR_EOF
true || echo 'restore of calc-comb.el failed'
fi
echo 'End of  part 10'
echo 'File calc-comb.el is continued in part 11'
echo 11 > _shar_seq_.tmp
exit 0
exit 0 # Just in case...
--
Kent Landfield                   INTERNET: [e-mail]
Sterling Software, IMD           UUCP:     uunet!sparky!kent
Phone:    [phone]         FAX:      [phone]
Please send comp.sources.misc-related mail to [e-mail].